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We investigate the zero temperature phase diagram of a gas of bosonic atoms in one- and two- 
color standing-wave lattices in the framework of the Bose-Hubbard model. We first introduce some 
relevant physical quantities; superfluid fraction, condensate fraction, quasimomentum distribution, 
and matter-wave interference pattern. We then discuss the relationships between them on the formal 
level and show that the superfluid fraction, which is the relevant order parameter for the superfluid 
to Mott-insulator transition, cannot be probed directly via the matter wave interference patterns. 
The formal considerations are supported by exact numerical solutions of the Bose-Hubbard model 
for uniform one-dimensional systems. We then map out the phase diagram of bosons in non-uniform 
lattices. The emphasis is on optical two-color superlattices which exhibit a sinusoidal modulation of 
the well depth and can be easily realized experimentally. From the study of the superfluid fraction, 
the energy gap, and other quantities we identify new zero-temperature phases, including a localized 
and a quasi Bose-glass phase, and discuss prospects for their experimental observation. 



PACS numbers: 03.75.Kk, 03.75.Lm, 05.30.Jp, 73.43.Nq 



I. INTRODUCTION 

With the first experimental observation of the super- 
fluid to Mott-insulator transition in ultracold atomic 
Bose gases trapped in optical lattice potentials P, Q the 
unique prospects of this new class of systems became 
evident. Compared to experiments with Bose-Einstein 
condensates in shallow magnetic traps, the use of opti- 
cal standing-wave lattices gives access to a new regime 
which is dominated by strong correlations. The most 
striking ramiflcation of these correlations are quantum 
phase transitions which are driven by the interplay 
of the different contributions to the Hamiltonian of the 
system. 

The superfluid to Mott-insulator transition is governed 
by competition between kinetic energy and repulsive 
atom-atom interaction: in the presence of strong repul- 
sive interactions the ground state is almost a pure Fock 
state with a deflnite occupation number at the individual 
lattice sites. This is associated with a vanishing of the 
superfluid density which is the natural orderparameter 
for superfluid to insulator phase transitions [j, H ■ 

New quantum phase transitions can be expected if 
other terms in the Hamiltonian come into play. One in- 
teresting case is the presence of some additional spatial 
structure of the lattice potential. This might range from 
a regular modulation of the well depths to random dis- 
order. From several studies in the context of solid-state 
systems it is known that random disorder gives rise to 
new zero-temperature phases such as the Anderson lo- 
calized phase and the Bose-glass phase 0, Q . 

Cold atoms in optical lattices offer unparalleled possi- 
bilities to study these disorder-related phenomena Q . By 
a superposition of multiple laser beams it is possible to 
generate a huge variety of different lattice topologies in a 



perfectly controlled manner. At the same time Feshbach 
resonances can be utilized to modify the strengths of the 
interatomic interaction. These powerful methods of tun- 
ing the relevant parameters are supplemented by versatile 
techniques to probe the state of the many-boson system, 
e.g. by observing the matter- wave interference pattern 
after the atoms were released from the lattice 0, Q or by 
measuring the energy gap between the ground state and 
the first excited state |2| . Atoms in an optical lattice are 
therefore an ideal model system for the experimental in- 
vestigation of the fundamental questions associated with 
quantum phase transitions in lattice systems. 

The aim of this paper is twofold: Firstly, we want to 
review the important physical quantities which describe 
the properties of Bose gases in optical lattices and inves- 
tigate the relations between them, both, analytically and 
numerically. The emphasis in this part of the paper is on 
the superfluid fraction and its connection to directly ob- 
servable quantities. Secondly, we want to systematically 
explore the zero-temperature phase diagram for bosons 
in two-color superlattices. 

The theoretical framework for the description of 
strongly correlated Bose gases in lattices is the Bose- 
Hubbard model which is summarized in Sec. ITU 
We define the fundamental quantities, condensate frac- 
tion, quasimomentum distribution, and superfluid frac- 
tion, and discuss their relation on the formal level in Sec. 
IIHI The connection of these fundamental properties to 
the experimental observables, e.g., the matter- wave in- 
terference pattern, is investigated in Sec. IIVI The for- 
mal deflnition shows that the superfluid fraction depends 
crucially on the excitation spectrum The basic ex- 
perimental observable, i.e. the interference pattern, only 
probes the ground state of the system and cannot there- 
fore provide full information on the superfluid properties 
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of the system. Using the example of the superfluid to 
Mott-insulator transition we illustrate these differences in 
the exact numerical solution of the Bose-Hubbard model 
for a one-dimensional uniform lattice. Finally, in a Sec. 
IVI we extend our studies to non-uniform lattice poten- 
tials. The emphasis is on so-called two-color superlattices 
which exhibit a sinusoidal variation of the well depths 
0. We map out the zero-temperature phase diagram as 
function of the interaction strength and the amplitude 
of the modulation and compare the behavior of number 
fluctuations, condensate fraction, superfluid fraction and 
energy gap. This enables us to identify additional zero- 
temperature phases, i.e. a localized phase and a quasi 
Bose-glass phase, and to discuss possibilities for their ex- 
perimental detection. 

II. BOSE-HUBBARD MODEL 

Quantum phase transitions of Bose gases in optical 
lattices are associated with complex correlations in the 
many-body state which go far beyond a simple Gross- 
Pitaevskii-like description. A theoretical model which is 
capable of describing these dominating correlations is the 
Bose-Hubbard model. First introduced for model studies 
related to ^He liq uids in porous media and granular su- 
perconductors [lOllllI it was recently applied to ultracold 
atomic Bose gases in optical lattices 

In the following we will briefly review the assumptions 
which lead to the Bose-Hubbard Hamiltonian and sum- 
marize solution methods and the basic observables. Be- 
cause our emphasis is on conceptual and general aspects 
we restrict ourselves to one-dimensional lattices, the for- 
mal generalization to two or three dimensions is straight- 
forward. 



A. Bose-Hubbard Hamiltonian 

The basic assumption of the Bose-Hubbard model is 
that the lattice wells are sufficiently deep for the state 
of the system to be described using a basis of single- 
particle wavefunctions, localized at the individual lattice 
sites. Only these localized ground states are taken into 
account and excited vibrational states are neglected. In 
the language of band structure theory the model space 
of the Bose-Hubbard model comprises the lowest energy 
band only, all excited bands are excluded. The localized 
single-particles states are given by the Wannier functions 
for the lowest band. The energy gap between the ground 
state band and the first excited band has therefore to be 
sufficiently large for admixtures from the excited bands 
to be negligible. 

A natural way to characterize the many-boson states in 
the model space of the Bose-Hubbard model is the occu- 
pation number representation. Let's assume we are deal- 
ing with a system of N bosons in a lattice composed of / 
lattice sites. We introduce a set of occupation numbers 



{ni, n/} which specify the number of bosons in the lo- 
calized single-particle state at the individual lattice sites. 
The set of Fock states |ni, nj) for all possible compo- 
sitions of the occupation numbers under the constraint 
Si=i = ^ forms a complete basis the Bose-Hubbard 
model space. 

We can define the associated creation and annihilation 
operators, a| and a^, for a boson localized at the ith 
lattice site: 

a| \ni,...,ni, ...,ni) = \/n~+T \ni, ...,ni + l, ...,nj) 
ai |ni, ni, n/) = y/n^ |ni, - 1, n/) . 

The occupation number operator for the ith site is given 

by 

hi = aja^ . (2) 

Throughout this paper the indices i,j = 1, ...,/ label the 
lattice sites. All operator-valued quantities are marked 
by a hat. 

Using these creation and annihilation operators one 
can easily translate the many-body Hamiltonian of the 
system — consisting of the kinetic energy, the external lat- 
tice potential, and the two-body interaction — into its sec- 
ond quantized form. This procedure leads directly to the 
Bose-Hubbard Hamiltonian for the one-dimensional lat- 
tice system 

Ho = - J^ial^i&i + h.a.) 

^ e, hi — ^ h,(hj - 1) . 
1=1 1=1 

The first term describes the tunneling between adjacent 
lattice sites characterized by a strength parameter J. It is 
associated with the kinetic energy part of the first quan- 
tized Hamiltonian. In general we will use periodic bound- 
ary conditions, i.e., tunneling between the first and the 
last lattice site is included (the site index / -I- 1 which 
appears in the summation is implicitly replaced by 1). 

The second term in gives an on-site single-particle 
energy which originates from the external potentials and 
the on-site part of the kinetic energy. The on-site ener- 
gies ei are constant for a translationally invariant lattice 
and might be set to zero. However, in the presence of 
an additional parabolic trapping potential or for irreg- 
ular lattices the vary with the site index i. Finally, 
the third term comprises the on-site two-body interac- 
tion characterized by the interaction strength V. 

The parameters J, e^, and V are given by the matrix 
elements of the terms of the first-quantized Hamiltonian 
the in the localized single particle states 0,[l2- In the fol- 
lowing we will use the tunneling parameter J as the unit 
of energy; the dimensionless ratios V/ J and e^/ J are var- 
ied to explore different regions of the zero-temperature 
phase diagram of the lattice system. As mentioned in 
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the introduction one can envisage being able to tune the 
parameters of the Hamiltonian JJ)) experimentally, either 
by changing the lattice geometry or by exploiting a Fes- 
hbach resonance. 

Besides the fundamental restriction of the Bose- 
Hubbard model to the ground state band some additional 
assumptions are used to construct the Hamiltonian 
(i) Only tunneling between adjacent sites is included, 
long range hopping over several lattice sites is neglected. 
(a) Only short-range on-site interactions are included; 
long-range interactions as they can occur for example in 
the presence of dipole-dipole forces are not considered 
here, (m) The tunneling strength J is assumed to be in- 
dependent of the site index. This is obviously the case in 
regular lattice potentials, but it is an approximation for 
irregular lattices (13j . If required all these simplifications 
can be easily abandoned. 



B. Exact solution 

The most direct way to obtain the ground state and 
the excitation spectrum of a zero-temperature Bose gas 
in one-dimensional lattices is through the exact numer- 
ical diagonalization of the Bose-Hubbard Hamiltonian. 
Within the occupation number basis the construction of 
the matrix representation of the Bose-Hubbard Hamilto- 
nian is trivial. The only substantial hurdle is the dimen- 
sion of the number basis. For fixed number of sites / and 
number of particles N the dimension of the number basis 
is given by 



D 



{N + i-iy. 
m{i-iy. 



(4) 



It grows dramatically with increasing system size: For 
fixed filling factor N/I = 1 the basis dimension for / = 8 
is D ^ 6435, for /= 10 it grows to D = 92378, and for 
/ = 12 it reaches D — 1352078. The Hamilton matrix 
is, however, extremely sparse because only the tunneling 
term in Q generates off-diagonal matrix elements. Effi- 
cient iterative Lanczos-type algorithms can therefore be 
utilized to determine the few lowest eigenvalues and the 
corresponding eigenvectors. This allows to treat systems 
with up to / = TV = 12 on a standard PC. 

As the result of the numerical solution of the eigenvalue 
problem we obtain the energy eigenvalues E^'^'^ and the 

corresponding eigenvectors Ca^ for a few states {v = 
0, 1, 2...). The eigenvectors provide the coefficients of the 
expansion of the eigenstates in the number basis 



D 



(5) 



Here the index a = 1, ...,D labels the different Fock 
states, i.e., the different sets {ni, ...,ni}a of occupation 
numbers. 

Besides the exact numerical diagonalization, which ex- 
plicitly yields the interesting state vectors, it is possible 



to compute selected observables using Monte-Carlo sim- 
ulations 01 which allow to treat systems with several 
hundred lattices sites 0, . Because of its simplicity 
and transparency we will, nevertheless, restrict ourselves 
to the exact diagonalization for relatively small lattices, 
which already exhibit the important universal properties 
of this class of systems. 

Several approximation methods have been developed 
to solve the Bose-Hubbard model. One of them is the 
so-called mean-field approximation where the state is as- 
sumed to be a direct product of independent single-site 
states. Thereby complex correlations between the indi- 
vidual sites that can be of particular importance for the 
description of quantum phase transition cannot be de- 
scribed. On top of this independent-site approximation 
the Gutzwillcr ansatz 17j can be used to perform varia- 
tional calculations of the ground state. 

Another approximation scheme uses the Bogoliubov 
approach. Here the annihilation operators a^ are replaced 
by a complex amplitude Zi plus a fluctuation operator Si. 
The fluctuation part is assumed to be small such that 
terms involving squares of the fluctuation operators can 
be neglected. The application and the limitations of this 
approach are discussed in PJ| ■ 



C. Simple observables 

From the ground state |^'o) and a few excited states 
obtained from the exact diagonalization of the Bose- 
Hubbard Hamiltonian one can extract several simple ob- 
servables by computing expectation values. 

The simplest observable is the mean occupation num- 
ber of the different sites 



n, = (*o|h,: 1^- 



(6) 



For a translationally invariant lattice fii will be the same 
for all lattice sites and equal to the filling factor N/ 1 inde- 
pendent of the parameters V and J. This, however, does 
not mean that the structure of the state does not change. 
A more telling quantity are the fluctuations around this 
mean occupation number given by 



;*o|nn*o>-(*o|n^|*o> 



(7) 



These number fluctuations provide direct information on 
the structure of the ground state (jSj). If l^'o) is a super- 
position of many Fock states ...,ni}a) then many 
different occupation numbers occur at the individual sites 
and the number fluctuations ai will be large. If, on the 
other hand, the ground state is described by a single Fock 
state, then the fluctuations will vanish. 

As a complementary measure the magnitude of the 
largest coefficient in the Fock state expansion (jsj of the 
ground state can be used 



^max — max(C„) 



(8) 
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FIG. 1: Number fluctuations (Ji (a) and largest coeflicient 
Cmax (b) for the ground state of the one-dimensional Bose- 
Hubbard model with 7 = A*' = 12 (solid), 10 (dashed), and 
8 (dotted) as function of the interaction strength V/ J. The 
vertical gray line marks the critical interaction strength for 
the Mott-insulator transition in an infinite system ,14. .19,1. 

Although this is not directly observable it yields sensitive 
information on the structure of the state. A small value 
of C^ax indicates that is a superposition of many 

number states; C^^x ~ 1 means that the ground state is 
a pure Fock state. 

To illustrate these quantities we solve the eigenvalue 
problem of the Bose- Hubbard Hamiltonian ||2J) for a reg- 
ular lattice potential (e^ ~ 0) with filling factor N/I = 1 
for different ratios V/J of the interaction strength and 
the tunneling coefficient. The results for the number fluc- 
tuations Ui and the largest coefficient C^^y. are depicted 
in Fig. n] The vertical gray line indicates the critical 
interaction strength (F/J)cr = 4.65 for the superfluid to 
Mott-insulator transition in an infinite one-dimensional 
lattice extracted from a Monte Carlo calculation nJI and 
from a strong coupling expansion 19]. It marks the re- 
gion where structural changes in the ground state are also 
expected in the finite size systems considered here. 

Our first important observation is that the number 
fiuctuations decrease rather slowly with increasing inter- 
action strength. Moreover, there appears no clear signa- 
ture of a phase transition in their variation. This is not 
an artifact of the limited lattice size, because the num- 
ber fluctuation exhibit no noticeable size dependence as 
the three curves for I = 12, 10, and 8 in Fig. ^a) illus- 
trate. Thus we conclude that the number fluctuations 
and related quantities do not reveal the quantum phase 
transition from a superfluid to a Mott-insulator state and 
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FIG. 2: Energy gap -Egap between ground state and first ex- 
cited state as function of the interaction strength V j J for com- 
mensurate filling N/I — 1 and different lattice sizes 7 = 12 
(solid line), 10 (dashed), 8 (dotted). 

that the insulating state still has significant number fiuc- 
tuations, i.e., is not a pure Fock state. 

The last point is confirmed by the behavior of the 
largest coefficient C^^x plotted in Fig. ^b). This quan- 
tity is very small in the insulating phase and starts to 
grow in the region of the phase transition. However, even 
for large values of the ratio V/ J the largest coefficient in 
(|SJ| remains smaller than 1 indicating that there is a dom- 
inant Fock state (the one with rii = N/I) but that other 
Fock states also contribute to the ground state. 

Another important observable is the energy gap, i.e., 
the difference between the energies of the ground state 
and the first excited state, 

i?g,p = £;W_ij(o). (9) 

This quantity is of particular importance as the excited 
states are crucial for the superfluid properties of the sys- 
tem. Figure [21 depicts the behavior of the energy gap 
as function of the interaction strength. For small V/J 
the energy gap is small and increases only slowly with 
increasing interaction strengths. The slope then starts 
to increase until, for V/J > 10, one observes a linear in- 
crease of the energy gap with a slope 1, i.e., changing the 
interaction strength by a certain amount will change 
the energy gap by the same amount. Qualitatively, this 
behavior can be explained by assuming that the excited 
state is a superposition of different one-particle one-hole 
excitations of the ground state. The predominant Fock 
states in the expansion of the excited state thus have 
one doubly occupied site which leads to the linear depen- 
dence of the energy gap Eg^p on V. These signatures of 
the Mott-insulator phase have been investigated experi- 
mentally 2] as well as theoretically [23, |^ by tilting the 
lattice potential. 

The behavior of the energy gap indicates that the 
character of the energy spectrum changes in the region 
V/ J ^ 5 to 8. As we will see in Sec. 1111 HI this change 
is accompanied by a vanishing of superfluidity, i.e. the 
superfluid to Mott-insulator transition. 
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III. CONDENSATE AND SUPERFLUIDITY 

Besides the simple quantities discussed so far there are 
two more complex physical quantities which reveal sub- 
stantial information on fundamental physical properties 
of the Bose gas in the lattice. These are the condensate 
fraction and the superfluid fraction. This section is de- 
voted to the definition and discussion of these non-trivial 
quantities. 

A. Condensate fraction 




Unlike the Gross-Pitaevskii description of weakly in- 
teracting Bose gases the exact solution of the Bose- 
Hubbard model does not assume a perfect Bose-Einstein 
condensate from the outset. Moreover, due to the general 
representation of the ground state in terms of Fock states 
the presence or absence of a Bose-Einstein condensate is 
not obvious. For example, the exact ground state © of 
a noninteracting system at zero temperature, where one 
expects a perfect Bose-Einstein condensate, is a superpo- 
sition of many Fock states with coefficients given by the 
multinomial distribution. 



1. Definition 

In order to define a Bose-Einstein condensate for the 
ground state |\l/o) of a general interacting many-boson 
system in a lattice we consider the one-body density ma- 
trix associated with |^'o) 

p!j^ = (*o|a]a,|*o). (10) 

The eigenvectors of the one-body density matrix describe 
the so-called natural orbitals and the eigenvalues the cor- 
responding occupation numbers. 

Following the formulation of Penrose and Onsager 
a Bose-Einstein condensate is present if one of the nat- 
ural orbitals is macroscopically occupied. Its occupation 
number is just the number of condensate particles A'c and 
the eigenvector constitutes the condensate wave function 

Xc.i- 

Strictly speaking, macroscopic occupation of an orbital 
means that the ratio of its occupation number and the 
total particle number, in the following called condensate 
fraction fc ~ Nc/N, remains finite in the the thermody- 
namic limit (A^, / oo, N/I = const.). This implies that 
one cannot rigorously determine the absence of a Bose- 
Einstein condensate in finite-size systems. From the nor- 
malization of the one-body density matrix, Tr plj-* = N, 
it follows immediately that for any system there is an 
eigenvalue larger or equal N/I. Thus the condensate frac- 
tion fc is always larger than 1/J and cannot vanish for 
finite-size systems. Despite of these formal complications 
the ratio of the largest eigenvalue of (|10|l and the total 
particle number provides important physical information 



FIG. 3: Condensate fraction fc as function of the interaction 
strength V/ J for commensurate filling N/I = 1 and diflferent 
lattice sizes 1—12 (solid line), 10 (dashed), 8 (dotted). The 
thin lines at the right border indicate the asymptotic value 

/c ^ l/I- 



even in finite lattices. For simplicity we will use the term 
condensate fraction for this quantity. 

The presence of a Bose-Einstein condensate implies off- 
diagonal long range order That is, the matrix ele- 
ments of the one-body density matrix do not go to zero 
far off the diagonal but remain finite: 

p\f^O for \i-j\^oo. (11) 

The prove of this connection is straightforward if one con- 
siders the spectral decomposition of the one-body density 
matrix which allows us to separate the condensate part 
from a residual positive semi-definite density matrix p'j' 

= Ncx:,ac,,+f^ ■ (12) 

It is important to realize that off-diagonal long range or- 
der is associated directly with the phenomenon of Bose- 
Einstein condensation whereas the connection with su- 
perfluidity is, at most, indirect. 

The dependence of the condensate fraction, i.e., the rel- 
ative occupation of the natural orbital with largest eigen- 
value, on the interaction strength V/ J in a, regular lattice 
is depicted in Fig. O In the noninteracting case V/J = 
all particles are in the condensate. As soon as a finite re- 
pulsive interaction is present {V/J > 0) the condensate 
fraction is reduced, i.e., the condensate is depleted and 
non-condensate states are populated. For large interac- 
tion strengths the condensate fraction converges to l/I. 
This is the reason for the pronounced size dependence of 
the condensate fraction in the regime of strong depletion. 

2. Quasimomentum Distribution 

The eigenvalues of the one-body density matrix H1U|) 
provide more information than just the condensate frac- 
tion. For a translationally invariant lattice one can show 



6 



that the natural orbitals satisfy the Bloch theorem and 
therefore are quasimomentum eigenstates j25| . Thus the 
occupation numbers of the natural orbitals correspond 
to occupation numbers hq of the Bloch states with dif- 
ferent quasimomenta q. Notice that for a finite lattice 
of length L the quasimomenta can assume only discrete 
values which are integer multiples of 27r/L. The Bloch 
function for quasimomentum q — describes the conden- 
sate state as defined above. 

An alternative way to obtain the quasimomentum dis- 
tribution for a regular lattice is to use the relation be- 
tween the localized Wannier functions and delocalized 
Bloch functions, which are the quasimomentum eigen- 
states. The Bloch functions tpq{x) of the lowest band can 
be decomposed in terms of Wannier functions w{x — ^i): 



w{x - f,) , 



(13) 



where is the coordinate of the ith lattice site. Using the 
fact that the creation operators aj of the Bose-Hubbard 
model create a boson in the Wannier state w{x — S^i) we 
can readily define creation operators cj, f or a boson in the 
Bloch state with quasimomentum q |l2| 
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FIG. 4: Occupation numbers fig of Bloch states with quasi- 
momentum q function of the interaction strength V/J for a 
lattice with / = A'^ = 12. The solid line shows the condensate 
q = state, the dashed and dotted lines correspond to the 
first and second excited quasimomentum state, respectively. 



B. Superfluid fraction 

One of the most interesting physical quantities in an 
interacting many-boson system is the superfluid density. 
It is the natural order parameter for the various super- 
fluid to insulator phase transitions, e.g., the superfluid to 
Mott-insulator transition observed experimentally [l|. 



^ l—l 



(14) 



The occupation numbers of the Bloch states with quasi- 
momentum q are thus given by 



1 I 

\q = (vfol etc, iM'o) = J E e'«(«'-«^) {^o\ a]a. 

I 



(15) 



where a — ^i+i — £,i is the lattice spacing. As expected, 
the quasimomentum distribution is related to the Fourier 
transform of the one-body density matrix. 

Figure ^ shows the occupation number for the three 
lowest quasimomentum states as function of the interac- 
tion strength V/J. The population for the q = state 
corresponds to the condensate fraction and is depleted 
as the interaction strength is increased. Simultaneously 
with this, the occupation numbers of Bloch states with 
higher quasimomenta, which vanish in the noninteracting 
system, increase successively. The repulsive interactions 
thus lead to a redistribution of the particles from the con- 
densate to states with higher quasimomenta. In the limit 
of strong interactions all quasimomentum states exhibit 
the same occupation, i.e., the band is filled uniformly. 



1. Definition 

The macroscopic, as well as the microscopic definition 
of superfiuidity, are connected to fiow properties. Macro- 
scopically, one usually employs the two-fiuid picture to 
describe systems which exhibit superfiuidity. The super- 
fluid and the normal fluid are distinguished through their 
behavior in the presence of moving boundaries. For ex- 
ample, in a narrow channel with moving side walls only 
the normal fluid is dragged along with the walls, whereas 
the superfluid stays at rest. In the rest frame of the walls 
only the superfluid moves. From the difference between 
the energies (in the rest frame of the walls) of a sys- 
tem with moving and one with stationary walls, one can 
directly extract the kinetic energy of the superfluid com- 
ponent and hence, via its velocity, the superfluid density 

A microscopic definition of superfluidity can be con- 
structed in analogy to this simple macroscopic picture. 
First of all we consider the generalized condensate wave 
function x(af) — e"'^^^|x(a?)| as defined in the previous 
section. A spatial variation of the phase 'd{x) of the con- 
densate wave function is connected to a velocity field 



(16) 



The flow described by this velocity held is non-dissipative 
and irrotational, properties which conform with the 



7 



macroscopic picture of superfluidity. This velocity is 
therefore identified with the flow field of a superfluid [23| • 
In principle Eq. H16|l provides the basic microscopic dis- 
tinction of the superfluid through its flow behavior. 

As in the macroscopic case we can develop an explicit 
expression for the superfluid fraction by considering the 
total energy of the system. For simplicity we restrict 
ourselves to a one-dimensional system. Assume that we 
impose a linear phase variation with a total phase twist Q 
over the length L of the system, i.e., = Qx/L. One 
way to accomplished this is by imposing twisted bound- 
ary conditions on the many-body wave function 



^{xi, ...,Xk + L, 



...,a;fc, ...,a;Ar) (17) 



for all k |29l l30l |31| . This imposed phase gradient gener- 
ates a flow with a velocity Vg — HQ / (mL) . The portion of 
the system which responds to the phase gradient by flow- 
ing with this velocity is called superfluid. The presence 
of the flow will increase the total energy of the system. 
As long as other excitations are absent, that is, as long 
as the imposed phase gradient is small, the increase in 
the total energy can be attributed solely to the kinetic 
energy of the superflow 



Ee-Eo = ii\/s«s 



(18) 



Here Eq is the ground state energy of the system without 
phase twist and Eq is the ground state energy obtained 
with twisted boundary conditions. The mass Ms corre- 
sponds to the total mass of the superfluid portion, by 
introducing the superfluid fraction fs is can be rewritten 
as Ms — mN fs- From Eq. I|18(l we can extract a simple 
expression for the superfluid fraction [s^ [s^ 



2mL^ E, 







Eo 



Q2 



for e < TT . 



(19) 



Hence the superfluid fraction is determined by the stiff- 
ness of the system under phase variations. We have to 
stress that the twist angle Q has to be sufflciently small 
to avoid effects other than the collective flow of the su- 
perfluid component. We will return to this point at the 
end of this section. 

The definition of the superfluid fraction through the 
energy change under imposed phase variations corre- 
sponds to the so-called helicity modulus introduced in 
Ref. [23. Moreover, it is equivalent to the concept of 
winding numbers which is used in path-integral Monte 
Carlo methods 0| . All these measures predict the non- 
interacting Bose gas to be a perfect superfluid. This is 
in contradiction to the Landau picture of superfluidity 
which is based on the dispersion relation of elementary 
excitations which defines a critical velocity. In the Lan- 
dau picture the ideal Bose gas is not considered a super- 
fluid since its critical velocity is zero. A possible connec- 
tion between these two pictures of superfluidity is still a 
matter of debate 0|. We adopt a pragmatic view and 
use the definition H19|l keeping in mind that it does not 



tell anything about the dynamical stability of the super- 
fluid flow at finite velocities. 

A second comment is appropriate concerning the re- 
lation between superfluidity and Bose-Einstein conden- 
sation. First of all one should realize that although the 
phase of the condensate wave function determines the su- 
perfluid velocity the superfluid and the condensate 
fraction are not the same. A famous example is liquid 
^He at zero temperature which is 100% superfluid but 
only 10% of the atoms are in the condensate j^^. From 
the formal point of view we used the condensate wave 
function only to motivate the connection between twisted 
boundary conditions (|17|l and the superfluid flow. On the 
level of the resulting Eqs. (|17|l and (|19|l the condensate 
does not appear. In general it is not strictly evident that 
superfluidity presupposes the presence of a condensate 



The relation H19|l for the continuous system can be eas- 
ily transfered to the Bosc-Hubbard model. Basically one 
replaces the length L of the continuous system by the 
number of lattice sites / and the prcfactor of the kinetic 
energy li} / {2m) by the tunneling strength J of the Bose- 
Hubbard Hamiltonian. The superfluid fraction for the 
discrete lattice thus reads |5| 



fs 



P_ Ee 
JN 



Eo 



Q2 



for e < vr 



(20) 



In order to compute the ground state energy of a 
Bose-Hubbard system with imposed phase gradient, the 
use of twisted boundary conditions is impracticable. 
By means of the local unitary transformation Ue = 
exp(i8 Xk/L) we can map the phase twist from the 
state onto the Hamiltonian. The eigenvalues of the re- 
sulting twisted Hamiltonian for periodic boundary con- 
ditions are identical to those of the original Hamiltonian 
with twisted boundary conditions. The twisted Bose- 
Hubbard Hamiltonian reads 



Hf- 



J^(e-'«/^a|+ia,-Kh.a.) 



1=1 



n, -l--^^ni(ni-l) 



(21) 



i=l 



Compared to the original Hamiltonian © only the off- 
diagonal hopping term is modified by an additional factor 
gTi©/^^ the so-called Peierls phase factor 1^3, |^. Phys- 
ically this term ensures that any particle which tunnels 
to an adjacent site acquires the correct phase in order to 
establish the linear phase variation across the system. 

There are several ways to impose this phase factor ex- 
perimentally. One possibility is to add a linear external 
potential or a homogeneous electric field in the case of 
charged particles. Another possibility is a constant accel- 
eration of the lattice, which in the case of optical stand- 
ing wave lattices can be realized by ramping the detuning 
between the two counter-propagating beams. One might 
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FIG. 5: Difference Ee — Eo between the ground state energies 
with an imposed phase twist and without as function of the 
total twist angle O for a system with I — N — 8. The different 
curves correspond to V/J = (solid line), V/J = 3 (dashed), 
and V/J = 5 (dotted). The gray lines depict the energy 
difference for the first excited state of the twisted system. 



be able to utilize this possibility to probe superfluidity 
directly. 

A straightforward procedure to compute the superfluid 
fraction is now evident: We numerically solve the eigen- 
value problems of the twisted and the non- twisted Hamil- 
tonian. He and Hq, with periodic boundary conditions 
and insert the resulting ground state energies Eq and Eq 
into (PUjl . 

Before doing so, it is useful to study how the energy 
difference Eq — Eq entering into Eq. depends on 

the twist angle 8. We have already assumed that has 
to be sufficiently small to avoid excitations other than 
the collective flow of the superfluid component. Figure [S] 
depicts the result of a numerical calculation for Eq — Eq 
as function of Q for different values of V/ J . For the 
noninteracting system (black solid curve) the energy is 
proportional to 9^ up to large values of 9, i.e., the su- 
perfluid fraction is independent of the actual value of 9 
used in the calculation. However, one should note that 
at 9 = TT, which corresponds to antiperiodic boundary 
conditions, a level crossing between the ground state and 
the first excited state of the twisted system appears. The 
presence of interactions leads to a "level repulsion" be- 
tween these two states (see dashed and dotted curves in 
Fig. [SI. Hence the energy deviates significantly from 
the Eq — Eq (X behavior for 9 ~ tt. If one were to 
use Eq. (|2U|1 for 9 ~ tt one would underestimate the 
superfluid fraction significantly. A physically meaningful 
application of Eq. 1)2011 to determine the superfluid frac- 
tion therefore requires 9 <C tt. Formally one might even 
impose the limit 9 — s- 0. For all the following numerical 
calculations we will use 9 = 0.1. 




FIG. 6: (a) Superfluid fraction /s for one-dimensional Bose- 
Hubbard systems with I — N = 12 (solid), 10 (dashed), and 
8 (dotted) as function of the interaction strength V/J. (b) 
First and second order contribution to fs in the perturbative 
calculation. The vertical gray line marks the critical interac- 
tion strength {V/ J)ci ~ 4.65 for the Mott-insulator transition 
in an infinite system p^ll9| . 

2. Direct numerical calculation 

We now employ Eq. 1)20)1 to compute the superfluid 
fraction /g from the energies Eq and Eq obtained nu- 
merically. The results for systems with commensurate 
filling N/I — 1 and lattice sizes up to / = 12 are shown 
in Fig. infa). 

In the range 3 < V/J < 6 the superfluid fraction, i.e. 
the order parameter for the superfluid to Mott-insulator 
transition, declines rapidly with increasing interaction 
strength V/J. In these relatively small systems we al- 
ready observe a rapid transition from the superfluid phase 
at small V/J to an insulating phase with vanishing su- 
perfluid fraction at large V/J. The transition region is 
in good agreement with the critical point {V/ J)ct = 4.65 
for the inflnite lattice determined by Monte Carlo simu- 
lations Q or a strong coupling expansion [l9l |. 

The dependence of the superfluid fraction on the lat- 
tice size shown in Fig. EI a) is rather moderate. In the 
superfluid phase for up to y/J < 4 there is no notice- 
able size dependence. 9nly the shoulder of the curve ex- 
hibits a size dependence: With increasing / the decrease 
in fs steepens and the tail of the curve is compressed 
and shifted towards lower V/J. From this observation 
we conclude that systems with / « 10 already allow very 
good qualitative, and even semi-quantitative conclusions, 
about the superfluid behavior of larger systems. 
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The rapid decrease of the superfluid fraction is in 
stark contrast to the behavior of the number fluctuations 
shown in Fig. ^a). Evidently, the number fluctuations 
stifl exhibit values of up to ai « 0.3 in cases where the su- 
perfluid fraction has already vanished. This proves that 
the Mott-insulator state can still be a superposition of 
many Fock states and thus show considerable number 
fluctuations. Only in the limit V/ J oo when ai van- 
ishes it becomes a pure Fock state. 



3. Perturbative treatment 

Apart from the direct solution of the eigenvalue prob- 
lems of the twisted and the non-twisted Hamiltonian one 
can use a perturbative approach to determine the energy 
difference Eq — Eq. As the superfluid fraction H2U|) is de- 
fined only for very small twist angles O this perturbative 
treatment does not introduce any additional approxima- 
tion. All terms of the expansion that contribute in the 
limit 0^0 will be included. 

The derivation consists of two steps: We first expand 
the twisted Hamiltonian in a power series up to second 
order in O 



of the original Hamiltonian Hp, 



e ^ 62 . . 

He ^ Ho + y J - ^ T = Ho 



H, 



pert ; 



(22) 



where the current operator J and the kinetic energy op- 
erator T are given by: 



J — iJ^(a-_,_jai — h.a.) 

1=1 

T = -J^(al+ia,-Kh.a.) 



(23) 



We then compute the energy shift Eq — Eg caused by the 
perturbation Hpoit in H22|) via second order perturbation 
theory. The lowest order contributions to the energy shift 
are quadratic in 6. Inserting them into (|20|l leads to an 
expression for the superfluid fraction consisting of two 
terms 



fs — /s 



(1) 



(2) 



(24) 



The first order term is proportional to the kinetic energy 
expectation value 



f(i) - 



1 



^-0 T = 



(25) 



where Vfo is the ground state of the original Hamiltonian 
Ho and p^j'' the associated one-body density matrix intro- 
duced in Eq. (|10|l . The second order term involves ma- 
trix elements of the current operator between the ground 
state \'$q) and all the excited states j^"^) {v = 1,2,...) 



f(2) 



1 



J|*o) 



NJ 



(26) 



All higher-order terms of the perturbative expansion of 
the energy shift do not contribute to the superfluid frac- 
tion in the limit 8^0. Thus Eq. (|24|l is an exact 
expression for /g. 

We note that this derivation and the resulting expres- 
sion is closely related to the so-called Drude weight which 
characterizes the dc conductivity of charged fermionic 
systems |3^ . 

This formulation of the superfluid fraction provides a 
detailed insight into the mechanisms that govern the su- 
perfluid properties of the system. First of all it shows that 
superfluidity is not a property of the ground state, but 
rather the response of the system to an external pertur- 
bation. The second order term H26|l introduces an explicit 
dependence on the whole excitation spectrum. Since this 
second order contribution always lowers the total super- 
fluid fraction we can conclude that the coupling to excited 
states leads to a suppression of superfluidity. The first 
order term (|25|l . which is just the ground state expecta- 
tion value of the kinetic energy and can be expressed in 
terms of the one-body density matrix of the ground state 
alone, gives only an upper bound for fg [s^. 

The behavior of the first and second order contribu- 
tions to the superfluid fraction for a regular lattice are 
shown in Fig. Efb). The flrst order term H25|l . i.e. 
the rescaled kinetic energy expectation value, is approx- 
imately one for small V/J and decreases slowly with 
increasing interaction strength. Within the insulating 
phase, where the total superfluid fraction shown in Fig. 
|SJa) has vanished completely, the flrst order contribution 

/i^'' still has values of up to 0.3. 

The vanishing of the total superfluid fraction is caused 
by a characteristic behavior of the second order contribu- 
tion H2()l) : In the transition region fs shows a threshold- 
like increase from zero to its maximum value. For val- 
ues of V/J beyond this maximum the second order term 

quickly converges to the value of /i"'^' . Thus the vanish- 
ing of the superfluid fraction in insulating phases, e.g., 
the Mott-insulator, results from the cancellation between 
flrst and second order term which both have considerable 
size. 

This clearly demonstrates that the coupling to the ex- 
cited states and the structure of the excitation spectrum 
is essential for the superfluid properties of the system. 
The ground state alone, represented by the first order 
contribution fs^\ provides only limited information — 
namely an upper bound — on the superfluid fraction. 
Therefore, any observable which is sensitive to the ground 
state only, such as number fluctuations, condensate frac- 
tion, and interference pattern, cannot provide clear in- 
formation the superfluid properties of the system. 

The only other observable which is sensitive to the 
properties of the excitation spectrum is the energy gap 
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^-gap between the ground state and the first excited state. 
As we have shown in Sec. Ill CI the energy gap starts to 
grow linearly with V/ J in the Mott-insulator phase. This 
is another indicator for the change in the excitation spec- 
trum that also leads to the vanishing of the superfluid 
fraction. However, as Eq. H26(l illustrates, the depen- 
dence of the superfluid fraction on the excited states is 
more complex than just the energy gap. 



IV. EXPERIMENTAL OBSERVABLES 

The two quantities, condensate fraction and superfluid 
fraction, discussed in the preceeding section describe two 
important and fundamental properties of the system, 
however, they are not directly accessible to simple ex- 
periments. In this section we discuss two observables 
which are directly accessible to present day experiments 
and show their relation to the fundamental quantities 
discussed above. 



A. Matter- wave interference pattern 

The simplest experimental procedure to obtain infor- 
mation on the state of the Bose gas in the lattice is 
to release the atoms from the lattice and detect the 
matter- wave interference pattern after some time of flight 
0, 12, 13- The crucial question is, how much the inter- 
ference pattern can tell us about the superfluid or the 
condensed portion of the system. 

Neglecting interactions during the expansion of the gas 
we can describe the time evolution of the single-particle 
states, which are initially localized at the different lattice 
sites, by Gaussian wave-packets. The matter-wave inten- 
sity observed at a point y after ballistic expansion for a 
time r can be written as 0] 



X(y) = (^o|At(y)A(y)|vI/o) 
The amplitude operators are given by 

A(f) =^G»(y,r) a, , 

i=l 



(27) 



(28) 



where Gi {y, r) is the amplitude of a Gaussian wave- 
packet that originates from site i after free expansion 
for a time r. We are interested only in the generic struc- 
ture of the intensity pattern, i.e., the presence or the 
absence of interference peaks. Therefore, we can simplify 
the treatment significantly by dropping all the terms re- 
lated to the spatial envelope of the interference pattern. 
The Gaussian can thus be replaced by a phase factor and 
the amplitude operator reduces to 



1=1 



i=l 



,i S<p{y,T)i 



(29) 
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FIG. 7: Matter-wave intensity as function of the phase difler- 
ence Sep for a one- dimensional lattice system with / = A*' = 12 
and different interaction strengths V/J = (solid line), 5 
(dashed), and 8 (dotted). 



Here 0i (y , t) is the phase accumulated on the path from 
the ith site to the observation point y. In the far-field 
approximation this can be replaced by 6(f>{y^T)i, where 
S(j}{y, t) is the phase difference between paths originating 
from adjacent lattice sites. The matter-wave intensity 
(|27|l as function of 6(j) now reads 

W) = J E c''*(*-^') (vI/o|a]a,|*o) 



I 



(30) 



1 ^e'.0(.-.) ,(1) 



This result proves the direct relation between the inter- 
ference pattern and the quasimomentum distribution of 
the gas in the lattice given by Eq. l(TB|l . The occupa- 
tion number fiq of the quasimomentum q state is given 
by the intensity I{S(j}) for a phase difference Scj) — qa: 



hq — 2{S(I) — qa) 



(31) 



As a special case of this result, the condensate fraction, 
i.e. the occupation of the q = state, is proportional 
to the intensity for Scj) = 0. Of course, this result agrees 
with the naive view that the interference pattern sirngly 
gives the momentum distribution of the trapped gas {l3| . 

Furthermore Eq. H30II enables us to draw an important 
conclusion about the relation between the matter-wave 
interference pattern and the superfiuid properties of the 
system. The interference pattern depends exclusively on 
the one-body density matrix plj^' for the ground state. 
The superfiuid fraction, in contrast, depends crucially on 
the excited states of the system as we have discussed in 
detail in Sec. IIII B 31 Thus the interference pattern can- 
not give complete information on the superfluid proper- 
ties of the system. It does not probe the physics that is 
crucial for superfluidity. 

The matter-wave interference patterns which result 
from Eq. H30|) for three different interaction strength 
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V/J are depicted in Fig. [7| For the noninteracting sys- 
tem (solid curve) it exhibits sharp interference peaks at 
5(j} = Q and 27r which correspond to the principal peaks 
observed experimentally 0|. Because we discarded the 
terms related to the envelope, the intensity T{54>) is sim- 
ply a 27r-periodic function in 5(j). With increasing inter- 
action strength V/ J (dashed and dotted curves in Fig. CJ 
the amplitude of the principal peaks is reduced and an 
"incoherent background" emerges. Effectively this leads 
to a gradual broadening of the interference peaks which 
was observed experimentally 1]. 

The emergence of the background has a straightfor- 
ward interpretation in terms of the quasimomentum dis- 
tribution. The gray arrows in Fig. [3 mark the values 
of 5(j) which correspond to the discrete quasimomenta 
q allowed in the lattice. According to H31|) the inten- 
sity at these points just gives the occupation numbers 
for the different quasimomentum states. Thus the deple- 
tion of the principal peak with increasing V/J just indi- 
cates the depletion of the condensate. The emergence of 
a background shows that states with nonvanishing quasi- 
momenta are successively populated — as we have already 
discussed in Sec. IIII A 21 In the limit of strong interac- 
tions all quasimomentum states of the band are occupied 
uniformly and the corresponding interference pattern is 
perfectly flat. 

The fact that the vanishing of the superfluid fraction 
is not associated with the vanishing of the interference 
fringes is highlighted by these numerical results. The 
values of V/J used to compute the three intensity dis- 
tributions in Fig. correspond to superfluid fractions 
of /s = 1 (solid curve), /s « 0.5 (dashed), and /s « 
(dotted). In the latter case, interference fringes are still 
clearly visible although the system is a perfect insulator, 
i.e., the superfluid fraction is zero. Thus, we see that with 
increasing V/ J the superfluid component vanishes much 
earlier than the interference fringes. In other words, the 
Mott insulator phase still exhibits a degree of phase co- 
herence. 

As a simple quantity to characterize the interference 
pattern one might introduce the visibility of the fringes. 
The standard definition of the fringe visibility 

.1 , -^max -'-mm 

~ I +1 ■ ^ ' 

-'-max I -'-mm 

only relies on the maximum and minimum values of the 
intensity. For a regular lattice with even number of sites 
we immediately find Imax = 2^(0) and Imin ~ ^{'^)- On 
the basis of (|31|l we can identify the maximum (mini- 
mum) intensity with the largest (smallest) quasimomen- 
tum occupation number. This definition of the visibility 
therefore provides a direct measure for the inhomogene- 
ity of the quasimomentum distribution: V is 1 as long as 
there are unoccupied quasimomentum states and it goes 
to zero if the band is filled uniformly. 

The dependence of the visibility measure on the in- 
teraction strength is depicted in Fig. |H| The visibility 
decreases only very slowly with increasing V/J and can 




0.2- 

qI : L . 1 : 1 . 

5 10 15 20 

V/J 

FIG. 8: Visibility V of the interference fringes as function of 
V/J for lattices with commensurate filling N/I = 1 and 1=12 
(solid), 1=10 (dashed), and 1 = 8 (dotted) lattice sites. 

have values up to V ~ 0.7 in the insulating phase where 
fs ~ 0. This confirms the observation that superfluidity 
vanishes much faster than the interference fringes and 
thus the coherence of the system 0. The fact that the 
curves in Fig. [Sjare practically independent of the lattice 
size proves that this statement also holds for large lat- 
tices and is not an artifact of the relatively limited size. 
With regard to the quasimomentum distribution the slow 
decrease of V indicates that very large ratios V/ J are re- 
quired to achieve a perfectly uniform occupation of the 
band. 

One can define other measures for the fringe visibility 
which are more sensitive to particular changes in the in- 
terference patterns. However, the most important and 
decisive measure for the coherence properties of the sys- 
tem is the interference pattern itself. 

B. Structure factor 

In the presence of a periodic modulation of the lat- 
tice potential the atom density distribution will develop 
a spatial structure Q. Systems of this type will be 
discussed in detail in Sec. These periodic density 
modulations should be observable using Bragg diffrac- 
tion of light from the trapped atomic gas. Experiments of 
this kind have already been performed using laser-cooled 
atomic gases in optical lattices 36, 37]. 

The quantity directly accessible through Bragg diffrac- 
tion experiments is the structure factor S{k) which is 
given by [s^ 

1 ^ 

^(^) ^pT. e"=(«'-«^)(*o| |*o) . (33) 

This quantity is sensitive to density-density correlations 
and is therefore a suitable indicator for the presence of di- 
agonal long-range order. For a homogeneous occupation 
of all lattice sites, present in the regular lattices discussed 
so far, the structure factor exhibits distinct peaks at in- 
teger values of ka/{2n). The presence of any additional 
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periodic structure in the atomic density distribution will 
lead to the appearance of new peaks at noninteger val- 
ues of fca/(27r). These signal the presence of diagonal 
long-range order. From the position and the relative am- 
plitude of the additional peaks in the structure factor 
one can derive detailed information on the spatial den- 
sity distribution. We will discuss this in detail in Sec. 



V. TWO-COLOR SUPERLATTICES 

Having discussed the relevant observables on a formal 
level and compared their behavior for the superfluid to 
Mott-insulator transition in a regular lattice we now want 
to explore another aspect of the phase diagram of the zero 
temperature Bose gas. Optical lattices offer the unique 
possibility of modifying the topology of the lattice in a 
flexible and perfectly controlled manner. They thus open 
new possibilities for studying many-body physics in non- 
uniform lattice potentials. In contrast to most solid-state 
systems, where the precise topology of the irregular lat- 
tice is difficult to specify, the high degree of experimental 
control facilitates investigations into specific aspects of 
"disorder"-induced phenomena. 




FIG. 9: (a) Lattice potential generated by superposition of 
two optical standing waves (see text), (b) Corresponding dis- 
tribution of the on-site energies ei with modulation amplitude 
A (after subtraction of an over-all constant). 



tial also influences the tunneling coefficient J which, in 
general, varies for the different lattice sites. For the lat- 
tices considered here the influence of this variation can 
be neglected in good approximation [T^ . 



A. Optical two-color lattices 

The simplest way to experimentally generate a non- 
uniform lattice potential is through a superposition 
of two collinear optical standing waves with different 
wavelengths — a so called two-color lattice. The interfer- 
ence between these standing waves produces a superlat- 
tice with a sinusoidal modulation of the depth of the lat- 
tice wells. More complex topologies have already been re- 
alized experimentally by using a superposition of several 
laser beams with different wave vectors 'sq]. Truly ran- 
dom lattices can be generated by superimposing a speckle 
pattern . Thus optical lattices cover the whole range 
from spatially modulated lattice potentials to disordered 
systems. 

In the following we consider a superlattice generated by 
the interference of two standing waves with wavelengths 
Ai and A2 = |Ai. The intensity of the second laser field 
is assumed to be much smaller than the intensity of the 
first laser, in this example I2 = 0.04/i. This combination 
results in a superlattice consisting of identical unit cells 
each composed of 5 elementary lattice sites. Figure El^a) 
depicts the lattice potential C/iat(a;) for a certain value of 
the relative phase. 

In the language of the Bose-Hubbard model this spatial 
modulation of the well depths corresponds to a modula- 
tion of the on-site energies e^. We assume a sinusoidal 
variation of the Ci in the range [—A, 0] as depicted in 
Fig. Elb). The parameter A characterizes the ampli- 
tude of the modulation, for A = we recover the regular 
lattice potential. The modulation of the lattice poten- 



B. Localized phase 

As a first step we consider a noninteracting Bose gas in 
the superlattice depicted in Fig. |^ We solve the eigen- 
value problem for the Bose-Hubbard Hamiltonian includ- 
ing two unit cells with filling factor 1, i.e., / = = 10. 
For V/ J — we successively increase the amplitude A 
of the spatial modulation of the on-site energies. The re- 
sulting distributions of mean occupation numbers fii and 
number fluctuations ai for the individual sites are de- 
picted in Fig. IIUI For a regular system, without any mod- 
ulation of the Ci, the lattice is occupied homogeneously. 
With increasing modulation amplitude A/ J the particles 
preferentially occupy those lattice sites with lower as 
one would classically expect. The degree of localization 
is governed by the competition between the on-site ener- 
gies and the tunneling term of the Bose-Hubbard Hamil- 
tonian (PJ. For sufficiently large A/ J only the site with 
the lowest in each unit cell is populated. The num- 
ber fluctuations cji behave accordingly: with increasing 
A/ J the number fluctuations at the deepest lattice well 
grow whereas the fluctuations at all other wells slowly 
decrease. In the limit of large A/ J the ground state is a 
superposition of all Fock states which have nonzero oc- 
cupation numbers only at the deepest well of each unit 
cell. 

There are several possibilities for detecting this local- 
ized phase in an experiment. The simplest signature ap- 
pears in the interference pattern after the expansion of 
the gas. Figure ITTI shows the intensity distribution that 
results from Eq. H3U|) for three different values of A/ J. 
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FIG. 10: Mean occupation numbers fii (a) and number fluc- 
tuations (Ti (b) for the different lattice sites i = 1, 10 as 
function of the amplitude A/ J of the sinusoidal modulation 
for the noninteracting gas. 

In the absence of the spatial modulation (upper panel) 
we observe the interference pattern of a regular lattice 
discussed in Sec. IIV Al With increasing A/ J, that is in- 
creasing degree of localization, additional peaks emerge. 
First the peaks next to initial principal peaks appear, 
then the next-order peaks emerge. In the limit of per- 
fect localization we obtain a regular interference pattern 
with a peak separation corresponding to the inverse of the 
number of elementary sites in a unit cell. In this manner 
the relative strength of the emerging peaks provides a di- 
rect experimental measure for the degree of localization 
in these superlattices. 

Another experimental observable that is sensitive to 
the degree of localization is the structure factor S{ka) 
introduced in Sec. IIV Bl We will come back to this quan- 
tity in the next section. 



C. Quasi Bose-glass phase 

It is intuitively clear that repulsive interactions will 
have a strong influence on the localized phase: A con- 
centration of many particles at a single lattice site will 
generate a large interaction energy contribution. Thus 
the interaction enforces delocalization of the particles. 

To illustrate the competition between lattice modu- 
lation and repulsive two-body interaction Fig. ^] de- 
picts the change of the mean occupation numbers 
and the number fluctuations at with increasing inter- 
action strength V/J for a fixed modulation amplitude 
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FIG. 11: Interference pattern after release from the lattice 
for different values of the amplitude A/ J of the sinusoidal 
modulation of the on-site energies. 



A/J = 50. 

The perfectly localized configuration persists only at 
very small interaction strengths. Already for V/ J >2 it 
is energetically favorable to reduce the maximum occu- 
pation number and thus the interaction energy by redis- 
tributing particles to sites with higher on-site energies ei. 
This region of interaction-induced delocalization is some- 
times termed Anderson- glass phase 6, 13]. Initially the 
changes in the mean occupation numbers happen contin- 
uously. However, for V/J > 10 this behavior changes. 
As shown in Fig. I12r a) there are extended intervals of 
interaction strengths where the mean occupation num- 
bers are constant and approximately integer. Between 
these regions of "stable" configurations there are narrow 
windows in which some of the fii change by roughly ±1. 
In the present example one of these almost steplike rear- 
rangements happens at V^/ J « 30. This behavior of the 
mean occupation numbers is also reflected in the number 
fluctuations depicted in Fig. EJb). Within the stable 
regions all cr^ are rather small, i.e., the dominant con- 
tribution to the expansion ((SJ of the ground state comes 
from a single Fock state. The steplike rearrangements are 
accompanied by an increase of the number fluctuations 
of those sites that change their mean occupation number. 

It is the subtle interplay between interaction and lat- 
tice irregularity that drives these successive rearrange- 
ments between configurations with integer occupation 
numbers. The sequence of rearrangements depends cru- 
cially on the topology of the lattice, i.e., the set of on-site 
energies. If the distribution of the becomes more homo- 
geneous (meaning that the number of different values of 
the ei increases) then the number of possible rearrange- 
ments increases and the stable regions shrink. In the 
limit of a lattice with random disorder this leads to the 
so-called Bose-glass phase @, 01 • For the superlattices 
with a rather small number of different on-site energies 
we will use the term quasi Bose-glass phase. The study of 
the superfluid fraction in Sec. IV Dl shows that the quasi 
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FIG. 12: Mean occupation numbers (a) and the number 
fluctuations Oi (b) for the individual lattice sites i = 1, 10 
as function of the interaction strength V/ J for fixed modula- 
tion amplitude A/J = 50. 



Bose-glass phase is insulating and that the energy gap is 
small. 

If the interaction strength reaches a value F = A a 
final and very sharp rearrangement happens and the sys- 
tem assumes a homogeneous occupation of all sites. This 
is the transition to the Mott-insulator phase where de- 
spite of the presence of the irregular lattice the ground 
state is dominated by the Fock state with uniform occu- 
pation numbers ni = N / 1. This transition between quasi 
Bose-glass and Mott-insulator phase happens indepen- 
dently of the lattice topology if the interaction strength 
V exceeds the amplitude A of the lattice modulation. 

The sequence of rearrangements in the quasi Bose-glass 
phase can be directly observed in experiment. The ob- 
servable which is most sensitive to the changes in the 
spatial distribution of the atoms in the superlattice is 
the structure factor S{ka) introduced in Sec. II V Bl For a 
state with a homogeneous population of the lattice sites, 
e.g., in the Mott-insulator phase, S{ka) shows peaks at 
integer values of /ca/(27r) only. If mean-occupation num- 
bers, i.e. the diagonal elements of the one-body density 
matrix, exhibit some periodic long-range structure then 
new peaks will appear at noninteger values of fca/(27r) 
from which information on the long-range order can be 
extracted. 

The localized phase in a superlattice at small inter- 
action strengths is the most pronounced realization of 
diagonal long-range order: For the superlattice discussed 




ka/{27T) 

FIG. 13: Structure factor S{ka) for different values of the 
interaction strength V/J and fixed amplitude A/J = 50. 
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FIG. 14: Amplitude of the first peak of the structure factor 
S{ka) located around ka « 2n/5 as function of the interaction 
strength V/J for fixed A/J = 50. 

here each 5th lattice site is occupied and all others are 
practically empty. The corresponding structure factor 
S{ka) is shown in the upper panel of Fig. ^| Peaks of 
almost equal height emerge at fca/(27r) = iy/5 {v integer) 
indicating the perfect localization at a single site in each 
unit cell. With increasing interaction strength V/J the 
particles are gradually redistributed to other lattice sites. 
This leads to a suppression of the peaks at non-integer 
values of fca/(27r) as the two lower panels Fig. 1131 illus- 
trate. Thus the amplitude of the peaks in S{ka) provides 
detailed information on the complicated spatial structure 
within the quasi Bose-glass phase. 

To highlight this point further. Fig. 1141 shows the de- 
pendence of the amplitude of the first peak of S{ka) lo- 
cated around fca/(27r) w 1/5 as function of V/J for fixed 
amplitude A/J of the lattice modulation. The parame- 
ters correspond to those used in Fig. W]\ and cover the 
full range from the localized phase to the Mott-insulator. 
Obviously each of the rearrangements of the mean occu- 
pation numbers within the quasi Bose-glass phase leaves 
a distinct signature in the peak amplitude. The regions 
of stable configurations as well as the steplike rearrange- 
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FIG. 15: Contour plots of the largest coefficient Cmax (a), the 
maximum number fluctuations (Tmax (b), and the condensate 
fraction /c (c) in the V-A plane for the sinusoidal lattice with 
/ = iV = 10. 



ments that are visible in Fig. I12f a'l directly show up in 
S{ka). Thus the structure factor, which is experimen- 
tally accessible through Bragg scattering of light from 
the bosons in the lattice, provides a comprehensive in- 
sight into the spatial density structure and its complex 
dependency on the interaction strength and the modula- 
tion amplitude. 



D. Phase diagrams in the V-A plane 

We shall now illustrate the behavior of the superfluid 
fraction, condensate fraction, and energy gap using two- 



dimensional phase diagrams in the V-A plane. This rep- 
resentation exposes the complicated interplay between 
lattice modulation and two-body interaction which leads 
to the rich phase structure. 

The gross phase structure is already revealed in the 
behavior of the largest coefficient C^j^x ™ the expan- 
sion (|SJ) of the ground state which is depicted in Fig. 
I15f a). The dark shadings correspond to values of C^^^x 
close to 1, i.e. ground states which are almost pure Fock 
states. The prominent region of large C^^y. for values 
V > A corresponds to the Mott-insulator phase. The 
lobe-like structures at large modulation amplitudes A be- 
long to the quasi Bose-glass phase. They correspond to 
the regions of stable configurations discussed in connec- 
tion with Fig. 1121 in which the mean occupation numbers 
are constant and approximately integral. The valleys of 
small C^ax which separate the lobes are connected to the 
rearrangements between different stable configurations. 
As mentioned earlier the number and position of these 
lobes depends on the topology of the lattice potential. 

In close connection with the behavior of Cj'^jax ^-^^ the 
maximum number fluctuations (Tmax = max(cri) depicted 
in Fig. ICT b). Values of C^ax close to 1 imply that the 
number fluctuations at all lattice sites are small. Thus 
the Mott-insulator phase and the lobes of the quasi Bose- 
glass phase appear as regions of small number fluctua- 
tions (light shadings). In contrast, the rearrangement 
valleys in the quasi Bose-glass phase, which are associ- 
ated with small C^axi show increased number fluctua- 
tions. Number fluctuations of the order 1 emerge only 
within a narrow band at small values of V/ J. For large 
values of A/ J this corresponds to the localized phase dis- 
cussed in Sec. IV Bl 

The behavior of the condensate fraction /c is shown 
in Fig. llSf c). In line with variations of the maximum 
number fluctuations, fc assumes very small values in the 
Mott-insulator phase. Within the quasi Bose-glass it de- 
creases monotonically with increasing V/J, i.e., it does 
not show the lobe structure as the number fluctuations. 
Values of the condensate fraction above 0.6 result only 
in a small area at small V/J and A/ J and a very nar- 
row stripe atV/J< 2. Thus only the localized phase and 
the expected superfluid phase exhibit sizeable condensate 
fractions. 

So far, these quantities are determined by the prop- 
erties of the ground state of Bose gas. Following the 
discussion in Sec. IIIIBI thev therefore do not fully deter- 
mine the superfluid properties of the system. Figure [TBI 
shows the contour plots for two quantities which probe 
the excitation spectrum, these are the superfluid frac- 
tion and the energy gap. The plot for superfluid fraction 
fs in Fig. llGf a) pinpoints the boundaries of the super- 
fluid phase. Large values of the superfluid fraction (dark 
shadings) result only in a restricted parameter region at 
small values of V and A. Everywhere else — that is, in 
the Mott-insulator, the quasi Bose-glass, and the local- 
ized phase — the system is a perfect insulator. In compar- 
ing the boundaries of the superfluid phase defined in Fig. 
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FIG. 16: Contour plots of the superfluid fraction /s (a) and 
the energy gap -Egap (b) for the sinusoidal lattice with / = 
N = 10. 



llGf a) with the structure of the other observables, Uke 
number fluctuations and condensate fraction in Fig. 1151 
it becomes evident that none of the ground state quan- 
tities can identify the superfluid region. Not only that 
the magnitude of the superfluid fraction cannot be deter- 
mined, the whole structure of the plots is different. 

This also applies to the energy gap -Egap plotted in 
Fig. I16r b). Although the opening of a gap between the 
ground state and the first excited state is an indicator for 
the superfluid to Mott-insulator transition in the uniform 
lattice, there is no such correspondence for the irregular 
lattice. It is only the Mott-insulator phase that exhibits a 
large energy gap. The localized and the quasi Bose-glass 
phase, which are also insulating, have small energy gaps 
just like the superfluid phase. 

The energy gap of the quasi Bose-glass phase deserves 
closer inspection. Figure I16r b') reveals that there is a 
small but finite energy gap within the stable lobes of the 
quasi Bose-glass phase. The size of the gaps, as well 
as the position and number of the lobes, depend on the 
differences between the various on-site energies. If the 
size of the unit cells is increased or if one goes to large 
random lattices then the differences between the be- 
come smaller, the number of lobes increases, their size 
decreases, and the energy gap vanishes. Thus there is 
a continuous transition from the quasi Bose-glass to the 
strict Bose-glass phase with vanishing -Egap- 



Another interesting observation results from the com- 
parison of the superfluid fraction in Fig. llGf a') with the 
condensate fraction in Fig. llSf c). In an irregular lat- 
tices with weak interactions there is a region where the 
condensate fraction is larger than the superfluid fraction. 
An extreme case is the localized phase at V^/ J < 2, where 
/s ~ and /c « 1. This means that the localized bosons 
in the different unit cells are still phase coherent but the 
superfluid flow is inhibited by the lattice topology. A sim- 
ilar phenomenon was recently observed in a Monte Carlo 
study of a continuous Bose gas with random impurities 
|40l | which act like the irregular lattice potential. 

So far we have only considered the phase diagrams for 
a filling factor N/I = 1. The qualitative structure of the 
phase diagram does not change if we go to larger inte- 
ger filling factors. Quantitatively the superfiuid region 
expands towards larger values of the interaction strength 
V/ J and the modulation amplitude A/ J if the filling fac- 
tor is increased. For a filling factor N/I — 4, for exam- 
ple, the contour corresponding to a superfluid fraction 
U = 0.6 extends up to V/J « 15 and A/J « 40. The 
insulating phases are shifted accordingly. 



E. Non-commensurate filling factors 

It is well-known that the Mott-insulator phase in a uni- 
form lattice exists only for commensurate filling factors. 
If one adds a few particles on top of a commensurate con- 
figuration then these particles can be moved to any lattice 
site without changing the energy — even in the presence of 
strong repulsive interactions. In this way the additional 
particles gives rise to a finite superfluid fraction in the 
limit of large V/J. 

This intuitive picture is confirmed by the exact calcu- 
lations for the sinusoidal superlattice with / = 10 sites 
and = 12 particles. Contour plots for the superfluid 
fraction and the energy gap in the V-IS. plane are pre- 
sented in Fig. El For the uniform lattice (A/ J = 0) the 
superfluid fraction shown in panel (a) remains finite even 
at large V/J . 

If one starts to increase the amplitude A of the sinu- 
soidal modulation then the superfluid fraction decreases 
and eventually vanishes above A/J > 5. In this regime 
the growing lattice irregularity causes a pinning of the 
additional particles to the lattice sites with minimal on- 
site energies e^. The ground state has an increased oc- 
cupation number at the deepest lattice well of each unit 
cefl. 

The behavior of the energy gap shown in Fig. llTr b) 
confirms this picture. For strong repulsive interactions 
the energy gap vanishes at small amplitudes A/ J in stark 
contrast to the commensurate system. With increasing 
amplitude of the sinusoidal modulation the energy gap 
starts to grow linearly with A/J. In this simple case 
Egap corresponds to the energy required to remove the 
additional particle from the lattice site it is affixed to. 
The slope of the increase is proportional to the energy 
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FIG. 17: Contour plots of the superfluid fraction /s (a) and 
the energy gap -Egap (b) in the V-/S. plane for the sinusoidal 
lattice with 7 = 10 and non-commensurable particle number 
N = 12. 



difference between the two lowest e^. This behavior is 
equivalent to the structure of the lobes of the quasi Bose- 
glass phase. 

Thus, in general, the change from a commensurate to a 
non-commensurate filUng factor causes a vanishing of the 
Mott-insulator phase. It is replaced by a superfluid phase 
at small A/ J, where the additional particles constitute a 
superfluid on top of an Mott-insulating layer. At larger 
A/ J the additional particles are pinned to the deepest 
lattice wells, superfluidity ceases, and an additional lobe 
of the quasi Bose-glass phase emerges. All other regions 
of the phase diagram are qualitatively unaffected by the 
change from a commensurate to a non-commensurate fill- 
ing factor. 



VI. CONCLUSIONS 

We have shown that the superfluid fraction — the nat- 
ural order parameter for superfluid to insulator phase 
transitions — is determined by the response of the system 
to an external perturbation, i.e., a twist in the boundary 
conditions, or equivalently an additional phase factor, in 
the tunneling term of the Bose-Hubbard Hamiltonian. 
This means that the superfluid fraction necessarily de- 
pends on the full excitation spectrum of the lattice sys- 



tem. The formal manifestation of this fact is given in Eq. 
(|24|l , where the superfluid fraction is expressed as a sum 
of a flrst order term depending only of the ground state 
and a second order term which involves all excited states. 
The exact numerical solution of the Bose-Hubbard model 
for uniform one-dimensional lattices reveals that the sec- 
ond order term is crucial for the vanishing of the super- 
fluid fraction in the Mott-insulator phase. The first order 
contribution gives only an upper bound for the superfluid 
fraction and can have a rather large value in the insulat- 
ing phase. 



This allows important conclusions for the prospects 
of measuring the superfluid fraction and the details of 
superfluid to insulator phase transitions. Present ex- 
perimental observables are sensitive only to the ground 
state and therefore cannot probe the physics that governs 
the superfluid properties. The matter-wave interference 
pattern, for example, provides direct information on the 
quasimomentum distribution of the ground state but it 
does not measure the superfluid properties. In the case of 
the superfluid to Mott-insulator transition our numerical 
results reveal that the superfluid fraction vanishes much 
earlier than the interference fringes, i.e., an interference 
pattern is still visible in the Mott-insulator phase. This 
demonstrates that one has to clearly distinguish between 
superfluidity and coherence properties. 



In the second part of the paper we have employed these 
tools to explore the zero-temperature phase diagram of 
Bose gases in non-uniform lattice potentials. We have 
demonstrated that even a simple superlattice potential 
which results from the superposition of two standing- 
wave lattices with different wavelengths gives rise to a 
very rich phase diagram. As function of the interac- 
tion strengths V/ J and the amplitude A/ J of the spa- 
tial variation of the on-site energies two additional in- 
sulating phases can be identifled: a localized phase and 
a quasi Bose-glass phase. All insulating phases can be 
clearly distinguished through their characteristic signa- 
tures in matter- wave interference experiments or through 
the structure factor determined by Bragg diffraction of 
light. 



Our results further support the view that Bose gases in 
optical lattices are a versatile tool for studying quantum 
mechanical many-body phenomena in strongly correlated 
systems. They allow us to address fundamental questions 
such as the connection between superfluidity and Bose- 
Einstein condensation. Moreover, they enable controlled 
studies of disorder-induced phenomena which result from 
a subtle competition between kinetic energy, two-body 
interaction, and lattice topology. To this end two-color 
lattices with relatively small unit cells seem a promising 
starting point because they facilitate a clear experimental 
distinction between the various phases. 
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